#include "homomorphic_encoding.h"
#include "../evaluator_utils.h"

#include <utility>

constexpr uint64_t FOO_MAX_LEN = UINT64_MAX; // 定义最大长度常量
/*
 * @brief EncodingMatrix 的构造函数。
 *
 * 初始化 EncodingMatrix 对象，根据给定的秘密密钥和方案参数创建编码矩阵。
 *
 * @param secretKey 秘密密钥。
 * @param scheme 方案参数。
 */
EncodingMatrix::EncodingMatrix(SecretKey &secretKey, Scheme &scheme) : secretKey(secretKey), scheme(scheme)
{
    long slots = scheme.context.N / 2;
    vector<complex<double>> A = ComputeRoots(slots << 1);
    vector<uint32_t> rotGroup(slots);

    rotGroup[0] = 1;
    for (int i = 1; i < slots; i++) {
        rotGroup[i] = rotGroup[i - 1] * 5;
        rotGroup[i] &= (slots << 2) - 1;
    }

    int levelBudgetEnc = 3;
    int levelBudgetDec = 3;

    mParamsEnc = GetCollapsedFFTParams(slots, levelBudgetEnc, 0);
    mParamsDec = GetCollapsedFFTParams(slots, levelBudgetDec, 0);

    mU0hatTPreFFT = EvalCoeffsToSlotsPrecompute(A, rotGroup, false, scheme.context.p, scheme.context.L - 1);
    mU0PreFFT = EvalSlotsToCoeffsPrecompute(A, rotGroup, false, scheme.context.p, scheme.context.L - 1 - levelBudgetDec - 9);
}

/*
 * @brief ：计算单位根数组，用于傅里叶变换或其他编码操作。
 *
 * 该函数计算 `2 * N` 个单位根，这些根分布在复数平面单位圆上，用于快速傅里叶变换（FFT）等算法。
 * 每个根的角度按 `2πi / (2 * N)` 计算，其中 `i` 是根的索引。
 * 计算结果包含一个 `2 * N + 1` 长度的复数数组，数组的最后一个根与第一个根相同，形成闭环。
 *
 * @param N：计算单位根所基于的整数，决定根的数量（`2 * N` 个根）。
 * @param a：该参数目前未使用，可能是为未来扩展预留的。
 *
 * @return 返回一个 `std::vector<std::complex<double>>` 类型的数组，包含计算得到的 `2 * N` 个单位根。
 */
vector<complex<double>> EncodingMatrix::ComputeRoots(int N)
{
    int m = N << 1;
    vector<complex<double>> roots(m + 1);
    roots[0] = std::complex<double>(1, 0);
    double angle;
    for (int i = 1; i < m; i++) {
        angle = 2 * 3.14159265358979323846 * static_cast<double>(i) / static_cast<double>(m);
        roots[i] = std::complex<double>(std::cos(angle), std::sin(angle));
    }
    roots[m] = roots[0];
    return roots;
}

/*
 * @brief ：旋转输入数组，将其按给定索引进行旋转。
 *
 * 该函数通过根据指定的索引旋转输入数组 `a`，生成一个新的数组。旋转会通过从数组的指定位置移动元素并将其放置在数组的另一端来实现。
 * 如果旋转索引超出了数组大小范围，会调用 `ReduceRotation` 方法来处理。
 *
 * @param a：输入的复数数组，包含需要旋转的数据。
 * @param index：旋转的索引，决定数组中各元素的位置变化。
 * @param tr：一个布尔值参数，目前未使用，可能是为未来扩展保留的。
 *
 * @return 返回旋转后的复数数组，大小与输入数组相同。
 */
complex<double> *EncodingMatrix::Rotate(const vector<complex<double>> &a, int index)
{
    int32_t slots = a.size();
    if (slots == 0 || slots > FOO_MAX_LEN) {
        throw std::runtime_error("Invalid value for slots.");
    }
    auto result = new complex<double>[slots]();

    if (index < 0 || index > slots) {
        index = ReduceRotation(index, slots);
    }

    if (index == 0) {
        for (int i = 0; i < slots; i++) {
            result[i] = a[i];
        }
    } else {
        // two cases: i+index <= slots and i+index > slots
        for (int32_t i = 0; i < slots - index; i++) {
            result[i] = a[i + index];
        }
        for (int32_t i = slots - index; i < slots; i++) {
            result[i] = a[i + index - slots];
        }
    }

    return result;
}

/*
 * @brief ：选择层数以适应给定的预算。
 * 
 * 根据给定的 `logSlots`（槽的对数值）和 `budget`（预算），该函数计算如何将输入的槽数量分配到多个层级中，以便满足预算要求。
 * 返回的结果表示所需的层数、每层的行数和剩余层数。
 * 
 * @param logSlots：表示槽数的对数值，用于计算层数。
 * @param budget：指定每个级别的预算，决定每个层级的大小。
 * 
 * @return 返回一个包含 `layers`、`rows` 和 `rem` 的向量，分别表示层数、每层的行数和剩余层数。
 */
vector<uint32_t> EncodingMatrix::SelectLayers(uint32_t logSlots, uint32_t budget)
{
    if (budget == 0) {
        throw std::runtime_error("Division by zero is not allowed.");
    }
    uint32_t layers = std::ceil(static_cast<double>(logSlots) / budget);
    if (layers == 0) {
        throw std::runtime_error("Division by zero is not allowed.");
    }
    uint32_t rows = logSlots / layers;
    uint32_t rem = logSlots % layers;

    uint32_t dim = rows;
    if (rem != 0) {
        dim = rows + 1;
    }

    // the above choice ensures dim <= budget
    if (dim < budget) {
        layers -= 1;
        rows = logSlots / layers;
        rem = logSlots - rows * layers;
        dim = rows;

        if (rem != 0) {
            dim = rows + 1;
        }

        // the above choice endures dim >=budget
        while (dim != budget) {
            rows -= 1;
            rem = logSlots - rows * layers;
            dim = rows;
            if (rem != 0) {
                dim = rows + 1;
            }
        }
    }

    return { layers, rows, rem };
}

/*
 * @brief ：根据给定的旋转组和预算层级计算编码系数。
 * 
 * 该函数根据 `rotGroup` 和 `levelBudget` 参数计算编码系数，并在多个层级中进行合并。首先通过 `CoeffEncodingOneLevel` 获取一层的编码系数，
 * 然后根据预算层级将这些系数展开并进行旋转处理。此函数还处理剩余层数并确保符合预算要求。
 * 
 * @param pows：用于计算编码系数的复数数组，表示旋转操作的幂次。
 * @param rotGroup：旋转组，表示每个层级的旋转操作。
 * @param levelBudget：预算层级数，决定最终计算的层级数量。
 * @param flag_i：布尔标志，指示是否采用特定的 `i` 处理。
 * 
 * @return 返回一个三维复数数组，表示不同层级的编码系数。
 */
std::vector<std::vector<std::vector<std::complex<double>>>> EncodingMatrix::CoeffEncodingCollapse(
        vector<std::complex<double>> pows, vector<uint32_t> rotGroup, uint32_t levelBudget, bool flag_i)
{
    uint32_t slots = rotGroup.size();
    // Need to compute how many layers are collapsed in each of the level from the budget.
    // If there is no exact division between the maximum number of possible levels (log(slots)) and the
    // level budget, the last level will contain the remaining layers collapsed.
    int32_t layersCollapse;
    int32_t remCollapse;

    std::vector<uint32_t> dims = SelectLayers(std::log2(slots), levelBudget);
    layersCollapse = dims[0];
    remCollapse = dims[2];

    int32_t dimCollapse = int32_t(levelBudget);
    int32_t stop = 0;
    int32_t flagRem = 0;

    if (remCollapse == 0) {
        stop = -1;
        flagRem = 0;
    } else {
        stop = 0;
        flagRem = 1;
    }

    uint32_t numRotations = (1 << (layersCollapse + 1)) - 1;
    uint32_t numRotationsRem = (1 << (remCollapse + 1)) - 1;

    // Computing the coefficients for encoding for the given level budget
    std::vector<std::vector<std::complex<double>>> coeff1 = CoeffEncodingOneLevel(pows, rotGroup, flag_i);

    // Coeff stores the coefficients for the given budget of levels
    std::vector<std::vector<std::vector<std::complex<double>>>> coeff(dimCollapse);
    for (uint32_t i = 0; i < uint32_t(dimCollapse); i++) {
        if (flagRem) {
            if (i >= 1) {
                // after remainder
                coeff[i] = std::vector<std::vector<std::complex<double>>>(numRotations);
                for (uint32_t j = 0; j < numRotations; j++) {
                    coeff[i][j] = std::vector<std::complex<double>>(slots);
                }
            } else {
                // remainder corresponds to the first index in encoding and to the last one in decoding
                coeff[i] = std::vector<std::vector<std::complex<double>>>(numRotationsRem);
                for (uint32_t j = 0; j < numRotationsRem; j++) {
                    coeff[i][j] = std::vector<std::complex<double>>(slots);
                }
            }
        } else {
            std::vector<std::vector<std::complex<double>>> a(numRotations);
            coeff[i] = a;
            for (uint32_t j = 0; j < numRotations; j++) {
                coeff[i][j] = std::vector<std::complex<double>>(slots);
            }
        }
    }

    for (int32_t s = dimCollapse - 1; s > stop; s--) {
        int32_t top = int32_t(std::log2(slots)) - (dimCollapse - 1 - s) * layersCollapse - 1;

        for (int32_t l = 0; l < layersCollapse; l++) {
            if (l == 0) {
                coeff[s][0] = coeff1[top];
                coeff[s][1] = coeff1[top + std::log2(slots)];
                coeff[s][2] = coeff1[top + 2 * std::log2(slots)];
            } else {
                std::vector<std::vector<std::complex<double>>> temp = coeff[s];
                std::vector<std::vector<std::complex<double>>> zeros(numRotations,
                                                                     std::vector<std::complex<double>>(slots, 0.0));
                coeff[s] = zeros;
                uint32_t t = 0;

                for (int32_t u = 0; u < (1 << (l + 1)) - 1; u++) {
                    for (uint32_t k = 0; k < slots; k++) {
                        coeff[s][u + t][k] += coeff1[top - l][k] * temp[u][ReduceRotation(k - (1 << (top - l)), slots)];
                        coeff[s][u + t + 1][k] += coeff1[top - l + std::log2(slots)][k] * temp[u][k];
                        coeff[s][u + t + 2][k] += coeff1[top - l + 2 * std::log2(slots)][k] *
                                                  temp[u][ReduceRotation(k + (1 << (top - l)), slots)];
                    }
                    t += 1;
                }
            }
        }
    }

    if (flagRem) {
        int32_t s = 0;
        int32_t top = int32_t(std::log2(slots)) - (dimCollapse - 1 - s) * layersCollapse - 1;

        for (int32_t l = 0; l < remCollapse; l++) {
            if (l == 0) {
                coeff[s][0] = coeff1[top];
                coeff[s][1] = coeff1[top + std::log2(slots)];
                coeff[s][2] = coeff1[top + 2 * std::log2(slots)];
            } else {
                std::vector<std::vector<std::complex<double>>> temp = coeff[s];
                std::vector<std::vector<std::complex<double>>> zeros(numRotationsRem,
                                                                     std::vector<std::complex<double>>(slots, 0.0));
                coeff[s] = zeros;
                uint32_t t = 0;

                for (int32_t u = 0; u < (1 << (l + 1)) - 1; u++) {
                    for (uint32_t k = 0; k < slots; k++) {
                        coeff[s][u + t][k] += coeff1[top - l][k] * temp[u][ReduceRotation(k - (1 << (top - l)), slots)];
                        coeff[s][u + t + 1][k] += coeff1[top - l + std::log2(slots)][k] * temp[u][k];
                        coeff[s][u + t + 2][k] += coeff1[top - l + 2 * std::log2(slots)][k] *
                                                  temp[u][ReduceRotation(k + (1 << (top - l)), slots)];
                    }
                    t += 1;
                }
            }
        }
    }

    return coeff;
}

/*
 * @brief ：计算编码过程中的一层系数。
 * 
 * 该函数通过从输入的复数数组 `pows` 和旋转组 `rotGroup` 计算出每个层级的编码系数。函数通过递归的方式（使用二进制递减）逐步计算每一层的系数。
 * 该过程主要依赖快速傅里叶变换（FFT）的思想，将输入数组进行多次左移和右移操作。
 * 
 * @param pows：用于旋转的复数幂数组。
 * @param rotGroup：旋转组，表示每个层级的旋转。
 * @param flag_i：布尔标志，指示是否采用特殊处理（`i`）。
 * 
 * @return 返回一个二维复数数组，表示每一层的编码系数。
 */
std::vector<std::vector<std::complex<double>>> EncodingMatrix::CoeffEncodingOneLevel(vector<std::complex<double>> pows,
                                                                                     vector<uint32_t> rotGroup, bool flag_i)
{
    uint32_t dim = pows.size() - 1;
    uint32_t slots = rotGroup.size();

    complex<double> I(0, 1);

    // Each outer iteration from the FFT algorithm can be written a weighted sum of
    // three terms: the input shifted right by a power of two, the unshifted input,
    // and the input shifted left by a power of two. For each outer iteration
    // (log2(size) in total), the matrix coeff stores the coefficients in the
    // following order: the coefficients associated to the input shifted right,
    // the coefficients for the non-shifted input and the coefficients associated
    // to the input shifted left.
    std::vector<std::vector<std::complex<double>>> coeff(3 * std::log2(slots));

    for (uint32_t i = 0; i < 3 * std::log2(slots); i++) {
        coeff[i] = std::vector<std::complex<double>>(slots);
    }

    for (uint32_t m = slots; m > 1; m >>= 1) {
        uint32_t s = std::log2(m) - 1;

        for (uint32_t k = 0; k < slots; k += m) {
            uint32_t lenh = m >> 1;
            uint32_t lenq = m << 2;

            for (uint32_t j = 0; j < lenh; j++) {
                uint32_t jTwiddle = (lenq - (rotGroup[j] % lenq)) * (dim / lenq);

                if (flag_i && (m == 2)) {
                    std::complex<double> w = std::exp(-M_PI / 2 * I) * pows[jTwiddle];
                    coeff[s + std::log2(slots)][j + k] = std::exp(-M_PI / 2 * I);     // not shifted
                    coeff[s + 2 * std::log2(slots)][j + k] = std::exp(-M_PI / 2 * I); // shifted left
                    coeff[s + std::log2(slots)][j + k + lenh] = -w;                   // not shifted
                    coeff[s][j + k + lenh] = w;                                       // shifted right
                } else {
                    std::complex<double> w = pows[jTwiddle];
                    coeff[s + std::log2(slots)][j + k] = 1;         // not shifted
                    coeff[s + 2 * std::log2(slots)][j + k] = 1;     // shifted left
                    coeff[s + std::log2(slots)][j + k + lenh] = -w; // not shifted
                    coeff[s][j + k + lenh] = w;                     // shifted right
                }
            }
        }
    }

    return coeff;
}

/*
 * @brief ：进行一次多项式系数解码，从输入数组 `pows` 和旋转组 `rotGroup` 计算出一个解码矩阵
 *          用于解码操作。
 * 
 * @param pows： 输入数组，包含模数集合中的多项式幂次。
 * @param rotGroup：旋转组数组，定义了旋转模式。
 * @param flag_i：布尔值标志，指示是否使用虚数单位（用于某些旋转计算）。
 * 
 * @return 返回一个二维复数矩阵，其中包含解码所需的系数。
 */
std::vector<std::vector<std::complex<double>>> EncodingMatrix::CoeffDecodingOneLevel(vector<std::complex<double>> pows,
                                                                                     vector<uint32_t> rotGroup, bool flag_i)
{
    uint32_t dim = pows.size() - 1;
    uint32_t slots = rotGroup.size();
    complex<double> I(0, 1);

    // Each outer iteration from the FFT algorithm can be written a weighted sum of
    // three terms: the input shifted right by a power of two, the unshifted input,
    // and the input shifted left by a power of two. For each outer iteration
    // (log2(size) in total), the matrix coeff stores the coefficients in the
    // following order: the coefficients associated to the input shifted right,
    // the coefficients for the non-shifted input and the coefficients associated
    // to the input shifted left.
    std::vector<std::vector<std::complex<double>>> coeff(3 * std::log2(slots));

    for (uint32_t i = 0; i < 3 * std::log2(slots); i++) {
        coeff[i] = std::vector<std::complex<double>>(slots);
    }

    for (uint32_t m = 2; m <= slots; m <<= 1) {
        uint32_t s = std::log2(m) - 1;

        for (uint32_t k = 0; k < slots; k += m) {
            uint32_t lenh = m >> 1;
            uint32_t lenq = m << 2;

            for (uint32_t j = 0; j < lenh; j++) {
                uint32_t jTwiddle = (rotGroup[j] % lenq) * (dim / lenq);
                if (flag_i && (m == 2)) {
                    std::complex<double> w = std::exp(M_PI / 2 * I) * pows[jTwiddle];
                    coeff[s + std::log2(slots)][j + k] = std::exp(M_PI / 2 * I); // not shifted
                    coeff[s + 2 * std::log2(slots)][j + k] = w;                  // shifted left
                    coeff[s + std::log2(slots)][j + k + lenh] = -w;              // not shifted
                    coeff[s][j + k + lenh] = std::exp(M_PI / 2 * I);             // shifted right
                } else {
                    std::complex<double> w = pows[jTwiddle];
                    coeff[s + std::log2(slots)][j + k] = 1;         // not shifted
                    coeff[s + 2 * std::log2(slots)][j + k] = w;     // shifted left
                    coeff[s + std::log2(slots)][j + k + lenh] = -w; // not shifted
                    coeff[s][j + k + lenh] = 1;                     // shifted right
                }
            }
        }
    }

    return coeff;
}

/*
 * @brief ：将多项式系数从一层展开到多个层次，并根据给定的层次预算压缩层级。
 * 
 * @param pows： 输入数组，包含多项式的幂次。
 * @param rotGroup：旋转组数组，定义了旋转方式。
 * @param levelBudget：层级预算，决定了将多少层进行压缩。
 * @param flag_i：布尔值，指示是否使用虚数单位（用于旋转计算）。
 * 
 * @return 返回一个三维复数矩阵，表示解码的系数。
 */
std::vector<std::vector<std::vector<std::complex<double>>>> EncodingMatrix::CoeffDecodingCollapse(
        std::vector<std::complex<double>> pows, std::vector<uint32_t> rotGroup, uint32_t levelBudget, bool flag_i)
{
    uint32_t slots = rotGroup.size();
    // Need to compute how many layers are collapsed in each of the level from the budget.
    // If there is no exact division between the maximum number of possible levels (log(slots)) and the
    // level budget, the last level will contain the remaining layers collapsed.
    int32_t layersCollapse;
    int32_t rowsCollapse;
    int32_t remCollapse;

    std::vector<uint32_t> dims = SelectLayers(std::log2(slots), levelBudget);
    layersCollapse = dims[0];
    rowsCollapse = dims[1];
    remCollapse = dims[2];

    int32_t dimCollapse = int32_t(levelBudget);
    int32_t flagRem = 0;

    if (remCollapse == 0) {
        flagRem = 0;
    } else {
        flagRem = 1;
    }

    uint32_t numRotations = (1 << (layersCollapse + 1)) - 1;
    uint32_t numRotationsRem = (1 << (remCollapse + 1)) - 1;

    // Computing the coefficients for decoding for the given level budget
    std::vector<std::vector<std::complex<double>>> coeff1 = CoeffDecodingOneLevel(pows, rotGroup, flag_i);

    // Coeff stores the coefficients for the given budget of levels
    std::vector<std::vector<std::vector<std::complex<double>>>> coeff(dimCollapse);

    for (uint32_t i = 0; i < uint32_t(dimCollapse); i++) {
        if (flagRem) {
            if (i < levelBudget - 1) {
                // before remainder
                coeff[i] = std::vector<std::vector<std::complex<double>>>(numRotations);

                for (uint32_t j = 0; j < numRotations; j++) {
                    coeff[i][j] = std::vector<std::complex<double>>(slots);
                }
            } else {
                // remainder corresponds to the first index in encoding and to the last one in decoding
                coeff[i] = std::vector<std::vector<std::complex<double>>>(numRotationsRem);

                for (uint32_t j = 0; j < numRotationsRem; j++) {
                    coeff[i][j] = std::vector<std::complex<double>>(slots);
                }
            }
        } else {
            coeff[i] = std::vector<std::vector<std::complex<double>>>(numRotations);

            for (uint32_t j = 0; j < numRotations; j++) {
                coeff[i][j] = std::vector<std::complex<double>>(slots);
            }
        }
    }

    for (int32_t s = 0; s < rowsCollapse; s++) {
        for (int32_t l = 0; l < layersCollapse; l++) {
            if (l == 0) {
                coeff[s][0] = coeff1[s * layersCollapse];
                coeff[s][1] = coeff1[std::log2(slots) + s * layersCollapse];
                coeff[s][2] = coeff1[2 * std::log2(slots) + s * layersCollapse];
            } else {
                std::vector<std::vector<std::complex<double>>> temp = coeff[s];
                std::vector<std::vector<std::complex<double>>> zeros(numRotations,
                                                                     std::vector<std::complex<double>>(slots, 0.0));
                coeff[s] = zeros;

                for (uint32_t t = 0; t < 3; t++) {
                    for (int32_t u = 0; u < (1 << (l + 1)) - 1; u++) {
                        for (uint32_t k = 0; k < slots; k++) {
                            if (t == 0)
                                coeff[s][u][k] += coeff1[s * layersCollapse + l][k] * temp[u][k];
                            if (t == 1)
                                coeff[s][u + (1 << l)][k] +=
                                        coeff1[s * layersCollapse + l + std::log2(slots)][k] * temp[u][k];
                            if (t == 2)
                                coeff[s][u + (1 << (l + 1))][k] +=
                                        coeff1[s * layersCollapse + l + 2 * std::log2(slots)][k] * temp[u][k];
                        }
                    }
                }
            }
        }
    }

    if (flagRem) {
        int32_t s = rowsCollapse;

        for (int32_t l = 0; l < remCollapse; l++) {
            if (l == 0) {
                coeff[s][0] = coeff1[s * layersCollapse];
                coeff[s][1] = coeff1[std::log2(slots) + s * layersCollapse];
                coeff[s][2] = coeff1[2 * std::log2(slots) + s * layersCollapse];
            } else {
                std::vector<std::vector<std::complex<double>>> temp = coeff[s];
                std::vector<std::vector<std::complex<double>>> zeros(numRotationsRem,
                                                                     std::vector<std::complex<double>>(slots, 0.0));
                coeff[s] = zeros;

                for (uint32_t t = 0; t < 3; t++) {
                    for (int32_t u = 0; u < (1 << (l + 1)) - 1; u++) {
                        for (uint32_t k = 0; k < slots; k++) {
                            if (t == 0)
                                coeff[s][u][k] += coeff1[s * layersCollapse + l][k] * temp[u][k];
                            if (t == 1)
                                coeff[s][u + (1 << l)][k] +=
                                        coeff1[s * layersCollapse + l + std::log2(slots)][k] * temp[u][k];
                            if (t == 2)
                                coeff[s][u + (1 << (l + 1))][k] +=
                                        coeff1[s * layersCollapse + l + 2 * std::log2(slots)][k] * temp[u][k];
                        }
                    }
                }
            }
        }
    }

    return coeff;
}

/*
 * @brief ：预计算编码系数矩阵并将其转换为槽位格式
 *          本函数主要用于根据给定的多项式系数、旋转组以及其它参数，通过一系列计算生成与槽位对应的明文向量。
 *          该过程通常涉及旋转、缩放、编码等操作。
 * 
 * @param A：包含多项式系数的复数数组。
 * @param rotGroup：旋转组数组，定义了旋转模式。
 * @param flag_i：布尔值，指示是否使用虚数单位（用于某些旋转计算）。
 * @param scale：缩放因子，用于调整最终结果的尺度。
 * @param L：控制编码过程中层数的参数，通常表示当前层次。
 * 
 * @return 返回一个二维向量 `result`，其中每个元素是一个 `Plaintext` 对象，表示对应的明文。
 */
vector<vector<Plaintext>> EncodingMatrix::EvalCoeffsToSlotsPrecompute(const vector<complex<double>> &A,
                                                                      const vector<uint32_t> &rotGroup, bool flag_i, double scale, uint32_t L)
{
    double Importscale = 0.0357143;

    int BASE_NUM_LEVELS_TO_DROP = 1;

    int levelBudget = mParamsEnc[0];
    int layersCollapse = mParamsEnc[1];
    int remCollapse = mParamsEnc[2];
    int numRotations = mParamsEnc[3];
    int b = mParamsEnc[4];
    int g = mParamsEnc[5];
    int numRotationsRem = mParamsEnc[6];
    int bRem = mParamsEnc[7];
    int gRem = mParamsEnc[8];

    int stop = -1;
    int flagRem = 0;

    long N = scheme.context.N;
    long M = scheme.context.M;
    long K = scheme.context.K;

    long slots = N / 2;

    if (remCollapse != 0) {
        stop = 0;
        flagRem = 1;
    }

    scale = scheme.context.p;
    // result is the rotated plaintext version of the coefficients
    vector<vector<Plaintext>> result(levelBudget);
    for (uint32_t i = 0; i < uint32_t(levelBudget); i++) {
        vector<vector<complex<double>>> a1;
        if (flagRem == 1 && i == 0) {
            // remainder corresponds to index 0 in encoding and to last index in decoding
            //            result[i] = vector<Plaintext>(numRotationsRem);
            for (int j = 0; j < numRotationsRem; j++) {
                vector<complex<double>> a2(slots, 0);

                auto *resultMid = new uint64_t[(L + stop - levelBudget + 1 + K) << scheme.context.logN]();
                Plaintext plain(resultMid, N, slots, L + stop - levelBudget + 1 + K);
                result[i].push_back(plain);
                a1.push_back(a2);
            }
        } else {
            //            result[i] = vector<Plaintext>(numRotations);
            for (int j = 0; j < numRotations; j++) {
                vector<complex<double>> a2(slots, 0);

                auto *resultMid = new uint64_t[(L + i - levelBudget + 1 + K) << scheme.context.logN]();
                Plaintext plain(resultMid, N, slots, L + i - levelBudget + 1 + K);
                result[i].push_back(plain);
                a1.push_back(a2);
            }
        }
        cts_yuan.push_back(a1);
    }

    if (slots == M / 4) {
        // ------------------------------------------------------------------------------
        // fully-packed mode
        // ------------------------------------------------------------------------------

        auto coeff = CoeffEncodingCollapse(A, rotGroup, levelBudget, flag_i);
        //        int sb = 0;
        for (int32_t s = levelBudget - 1; s > stop; s--) {
            double scaleReal = scale;
            for (int scaleI = 0; scaleI < scheme.context.L - L + levelBudget - s - 1; scaleI++) {
                scaleReal = scaleReal * scaleReal / scheme.context.qVec[scheme.context.L - scaleI - 1];
            }
            for (int32_t i = 0; i < b; i++) {
                for (int32_t j = 0; j < g; j++) {
                    if (g * i + j != int32_t(numRotations)) {
                        uint32_t rot =
                                ReduceRotation(-g * i * (1 << ((s - flagRem) * layersCollapse + remCollapse)), slots);
                        //                        if ((flagRem == 0) && (s == stop + 1)) {
                        //                            // do the scaling only at the last set of coefficients
                        //                            for (uint32_t k = 0; k < slots; k++) {
                        //                                coeff[s][g * i + j][k] *= scale;
                        //                            }
                        //                        }

                        auto rotateTemp = Rotate(coeff[s][g * i + j], rot);

                        for (int iVals = 0; iVals < slots; iVals++) {
                            cts_yuan[s][g * i + j][iVals] = rotateTemp[iVals];
                        }
                        scheme.context.EncodePQ(result[s][g * i + j].mx, rotateTemp, scaleReal, slots,
                                                L + s - levelBudget + 1, K);
                        //                        delete vals;
                    }
                }
            }
        }

        if (flagRem) {
            double scaleReal = scale;
            for (int scaleI = 0; scaleI < scheme.context.L - L + levelBudget - stop - 1; scaleI++) {
                scaleReal = scaleReal * scaleReal / scheme.context.qVec[scheme.context.L - scaleI - 1];
            }
            for (int32_t i = 0; i < bRem; i++) {
                for (int32_t j = 0; j < gRem; j++) {
                    if (gRem * i + j != int32_t(numRotationsRem)) {
                        uint32_t rot = ReduceRotation(-gRem * i, slots);
                        for (uint32_t k = 0; k < slots; k++) {
                            coeff[stop][gRem * i + j][k] *= Importscale;
                        }

                        auto rotateTemp = Rotate(coeff[stop][gRem * i + j], rot);
                        for (int iVals = 0; iVals < slots; iVals++) {
                            cts_yuan[stop][gRem * i + j][iVals] = rotateTemp[iVals];
                        }
                        scheme.context.EncodePQ(result[stop][gRem * i + j].mx, rotateTemp, scaleReal, slots,
                                                L + stop - levelBudget + 1, K);
                    }
                }
            }
        }
        cout << "Precompute Finished" << endl;
    } else {
        throw std::runtime_error("Precompute Failed !");
    }

    return result;
}

/*
 * @brief ：预计算从槽位到系数的转换
 *          本函数用于根据给定的多项式系数、旋转组及其他参数，将槽位数据转换为系数矩阵。
 *          这一过程涉及旋转、解码和缩放操作。
 * 
 * @param A：包含多项式系数的复数数组。
 * @param rotGroup：旋转组数组，定义了旋转模式。
 * @param flag_i：布尔值，指示是否使用虚数单位（用于某些旋转计算）。
 * @param scale：缩放因子，用于调整最终结果的尺度。
 * @param L：控制解码过程的层数参数，通常表示当前层次。
 * 
 * @return 返回一个二维向量 `result`，其中每个元素是一个 `Plaintext` 对象，表示经过解码后的系数。
 */
std::vector<std::vector<Plaintext>> EncodingMatrix::EvalSlotsToCoeffsPrecompute(const vector<std::complex<double>> &A,
                                                                                const vector<uint32_t> &rotGroup, bool flag_i, double scale, uint32_t L)
{
    int BASE_NUM_LEVELS_TO_DROP = 1;

    int levelBudget = mParamsDec[0];
    int layersCollapse = mParamsDec[1];
    int remCollapse = mParamsDec[2];
    int numRotations = mParamsDec[3];
    int b = mParamsDec[4];
    int g = mParamsDec[5];
    int numRotationsRem = mParamsDec[6];
    int bRem = mParamsDec[7];
    int gRem = mParamsDec[8];

    int flagRem = 0;
    long N = scheme.context.N;
    long M = scheme.context.M;
    long K = scheme.context.K;

    long slots = N / 2;

    if (remCollapse != 0) {
        flagRem = 1;
    }

    // result is the rotated plaintext version of coeff
    vector<vector<Plaintext>> result(levelBudget);
    for (uint32_t i = 0; i < uint32_t(levelBudget); i++) {
        if (flagRem == 1 && i == uint32_t(levelBudget - 1)) {
            // remainder corresponds to index 0 in encoding and to last index in decoding
            //            result[i] = vector<Plaintext>(numRotationsRem);
            for (int j = 0; j < numRotationsRem; j++) {
                auto *resultMid = new uint64_t[(L - i + K) << scheme.context.logN]();
                Plaintext plain(resultMid, N, slots, L - i + K);
                result[i].push_back(plain);
            }
        } else {
            //            result[i] = vector<Plaintext>(numRotations);
            for (int j = 0; j < numRotations; j++) {
                auto *resultMid = new uint64_t[(L - i + K) << scheme.context.logN]();
                Plaintext plain(resultMid, N, slots, L - i + K);
                result[i].push_back(plain);
            }
        }
    }
    auto coeff = CoeffDecodingCollapse(A, rotGroup, levelBudget, flag_i);

    for (int32_t s = 0; s < levelBudget - flagRem; s++) {
        for (int32_t i = 0; i < b; i++) {
            for (int32_t j = 0; j < g; j++) {
                if (g * i + j != int32_t(numRotations)) {
                    uint32_t rot = ReduceRotation(-g * i * (1 << (s * layersCollapse)), slots);
                    if ((flagRem == 0) && (s == levelBudget - flagRem - 1)) {
                        // do the scaling only at the last set of coefficients
                        //                            for (uint32_t k = 0; k < slots; k++) {
                        //                                coeff[s][g * i + j][k] *= scale;
                        //                            }
                    }

                    auto rotateTemp = Rotate(coeff[s][g * i + j], rot);
                    scheme.context.EncodePQ(result[s][g * i + j].mx, rotateTemp, scale, slots, L - s, K);
                }
            }
        }
    }
    if (flagRem) {
        int32_t s = levelBudget - flagRem;
        for (int32_t i = 0; i < bRem; i++) {
            for (int32_t j = 0; j < gRem; j++) {
                if (gRem * i + j != int32_t(numRotationsRem)) {
                    uint32_t rot = ReduceRotation(-gRem * i * (1 << (s * layersCollapse)), slots);
                    //                    for (uint32_t k = 0; k < slots; k++) {
                    //                        coeff[s][gRem * i + j][k] *= scale;
                    //                    }
                    auto rotateTemp = Rotate(coeff[s][gRem * i + j], rot);
                    scheme.context.EncodePQ(result[s][gRem * i + j].mx, rotateTemp, scale, slots, L - s, K);
                }
            }
        }
    }
    //    }

    return result;
}

/*
 * @brief ：对密文的ax项进行模降
 * 
 * @param cipher：输入密文对象。
 * 
 * @return 返回新的密文 `Ciphertext`。
 */
Ciphertext EncodingMatrix::KeySwitchDownFirst(Ciphertext &cipher)
{
    long K = scheme.context.K;
    uint64_t *ctax = new uint64_t[(cipher.l + K) * cipher.N]();
    std::copy(cipher.ax, cipher.ax + (cipher.l + K) * cipher.N, ctax);
    for (long index = 0; index < scheme.context.K; index++) {
        uint64_t *ai = ctax + (index << scheme.context.logN);
        scheme.context.PiINTTAndEqual(ai, index);
    }

    for (long index = K; index < cipher.l + K; index++) {
        uint64_t *ai = ctax + (index << scheme.context.logN);
        scheme.context.QiINTTAndEqual(ai, index - K);
    }

    uint64_t *ax = new uint64_t[cipher.l * cipher.N]();
    // ax = scheme.context.BaseConvertPQl2Ql(ctax, cipher.l);
    scheme.context.BaseConvertPQl2Ql(ax, ctax, cipher.l);
    scheme.context.NTTAndEqual(ax, cipher.l);
    return Ciphertext(ax, cipher.bx, cipher.N, cipher.slots, cipher.l, cipher.scale);
}

/*
 * @brief ：对密文进行模降
 * @param cipher：输入密文对象。
 * 
 * @return 返回新的密文 `Ciphertext`。
 */
Ciphertext EncodingMatrix::KeySwitchDown(Ciphertext &cipher)
{
    long K = scheme.context.K;
    uint64_t *ctax = new uint64_t[(cipher.l + K) * cipher.N]();
    uint64_t *ctbx = new uint64_t[(cipher.l + K) * cipher.N]();
    std::copy(cipher.ax, cipher.ax + (cipher.l + K) * cipher.N, ctax);
    std::copy(cipher.bx, cipher.bx + (cipher.l + K) * cipher.N, ctbx);
    for (long index = 0; index < scheme.context.K; index++) {
        uint64_t *ai = ctax + (index << scheme.context.logN);
        uint64_t *bi = ctbx + (index << scheme.context.logN);
        scheme.context.PiINTTAndEqual(ai, index);
        scheme.context.PiINTTAndEqual(bi, index);
    }

    for (long index = K; index < cipher.l + K; index++) {
        uint64_t *ai = ctax + (index << scheme.context.logN);
        uint64_t *bi = ctbx + (index << scheme.context.logN);
        scheme.context.QiINTTAndEqual(ai, index - K);
        scheme.context.QiINTTAndEqual(bi, index - K);
    }

    uint64_t *ax = new uint64_t[cipher.l * cipher.N]();
    uint64_t *bx = new uint64_t[cipher.l * cipher.N]();
    // ax = scheme.context.BaseConvertPQl2Ql(ctax, cipher.l);
    // bx = scheme.context.BaseConvertPQl2Ql(ctbx, cipher.l);
    scheme.context.BaseConvertPQl2Ql(ax, ctax, cipher.l);
    scheme.context.BaseConvertPQl2Ql(bx, ctbx, cipher.l);
    scheme.context.NTTAndEqual(ax, cipher.l);
    scheme.context.NTTAndEqual(bx, cipher.l);
    return Ciphertext(ax, bx, cipher.N, cipher.slots, cipher.l, cipher.scale);
}

/*
 * @brief ：反转一个区间内的元素
 *          该函数用于反转给定区间内的元素，使得序列的顺序发生反转。
 * 
 * @param start：区间的起始指针。
 * @param end：区间的结束指针。
 */
void reverseSection(uint64_t *start, uint64_t *end)
{
    while (start < end) {
        std::swap(*start, *end);
        start++;
        end--;
    }
}

/*
 * @brief ：交换密文的两个部分
 *          该函数通过反转三个部分的元素，完成密文部分的交换。先反转前 l*N 个元素，再反转最后 k*N 个元素，最后反转整个区间。
 * 
 * @param a：待交换的数据指针。
 * @param l：区间的长度 l。
 * @param k：要交换的第二部分长度 k。
 * @param N：数据的大小 N。
 */
void swapSections(uint64_t *a, int l, int k, int N)
{
    // 反转前 l*N 个元素
    reverseSection(a, a + l * N - 1);

    // 反转最后 k*N 个元素
    reverseSection(a + l * N, a + l * N + k * N - 1);

    // 反转整个区间，实现最终的交换
    reverseSection(a, a + l * N + k * N - 1);
}

/*
 * @brief ：将密文从 Q 模转换为 PQ 模

 * @param result：存储结果的 `Ciphertext` 对象。
 * @param cipher：输入密文对象。
 */
void EncodingMatrix::KeySwitchExt(Ciphertext &result, Ciphertext cipher)
{
    long l = cipher.l;
    long N = cipher.N;
    long K = scheme.context.K;
    long logN = scheme.context.logN;

    result = cipher;

    result.ax = new uint64_t[(cipher.l + K) * N]();
    result.bx = new uint64_t[(cipher.l + K) * N]();
    std::copy(cipher.ax, cipher.ax + cipher.l * N, result.ax + K * N);
    std::copy(cipher.bx, cipher.bx + cipher.l * N, result.bx + K * N);
}

/*
 * @brief ：预计算同态旋转操作
 *          该函数用于执行旋转操作的预计算过程。通过对密文的 `ax` 部分进行 INTT 变换，接着进行基转换和 TNTT 变换，计算旋转操作所需的结果。
 * 
 * @param cipher：输入密文对象。
 * 
 * @return 返回一个向量，其中包含经过旋转预处理的密文数据。
 */
vector<uint64_t *> EncodingMatrix::EvalFastRotationPrecompute(Ciphertext &cipher)
{
    vector<uint64_t *> result;
    long rotSlots;

    uint64_t *ax = new uint64_t[cipher.l << scheme.context.logN]();
    scheme.context.INTT(ax, cipher.ax, cipher.l);

    int beta = ceil(cipher.l / float(scheme.context.K));
    result = scheme.context.BaseConvertQl2T(ax, cipher.l);
    for (int j = 0; j < beta; j++) {
        scheme.context.TNTTAndEqual(result[j], scheme.context.T);
    }
    return result;
}

/*
 * @brief ：计算密文与切换密钥的内积结果
 
 * @param levelQ：密文层数。
 * @param aPolyTs：包含多项式的密文数据数组。
 * @param key：切换密钥。
 * 
 * @return 返回计算后的内积结果。
 */
uint64_t *EncodingMatrix::InnerProduct(int levelQ, vector<uint64_t *> aPolyTs, vector<vector<vector<uint64_t>>> key)
{
    int k = scheme.context.K;
    int blockLen = ceil((levelQ + k) / (float(scheme.context.gamma)));
    int beta = ceil(levelQ / (float(k)));
    uint64_t tN;
    uint64_t tjN;
    int slice;

    uint64_t *result = new uint64_t[(levelQ + k) << scheme.context.logN]();

    for (int i = 0; i < beta; i++) {
        for (int j = 0; j < blockLen; j++) {
            for (int t = 0; t < scheme.context.T; t++) {
                uint64x2_t tVec = vdupq_n_u64(scheme.context.tVec[t]);
                uint64x2_t twotVec = vdupq_n_u64(scheme.context.twotVec[t]);
                uint64x2_t trVec_high = vdupq_n_u64(scheme.context.trVecHigh[t]);
                uint64x2_t trVec_low = vdupq_n_u64(scheme.context.trVecLow[t]);
                tN = t * scheme.context.N;
                tjN = j * scheme.context.TN + tN;
                for (int n = 0; n < scheme.context.N; n += 4) {
                    uint64x2_t aVec1 = vld1q_u64(&aPolyTs[i][tN + n]);
                    uint64x2_t keyVec1 = vld1q_u64(&key[i][j][tN + n]);
                    uint64x2_t aVec2 = vld1q_u64(&aPolyTs[i][tN + n + 2]);
                    uint64x2_t keyVec2 = vld1q_u64(&key[i][j][tN + n + 2]);

                    uint64x2_t result1, result2;
                    MulModBarrettNeon(&result1, aVec1, keyVec1, tVec, twotVec, trVec_high, trVec_low);
                    MulModBarrettNeon(&result2, aVec2, keyVec2, tVec, twotVec, trVec_high, trVec_low);

                    if (i == 0) {
                        vst1q_u64(&scheme.context.polyTBootStrap[j][tN + n], result1);
                        vst1q_u64(&scheme.context.polyTBootStrap[j][tN + n + 2], result2);
                    } else {
                        uint64x2_t tmpVal1 = vld1q_u64(&scheme.context.polyTBootStrap[j][tN + n]);
                        uint64x2_t tmpVal2 = vld1q_u64(&scheme.context.polyTBootStrap[j][tN + n + 2]);

                        uint64x2_t sum1 = vaddq_u64(tmpVal1, result1);
                        uint64x2_t sum2 = vaddq_u64(tmpVal2, result2);

                        uint64x2_t cmp1 = vcgtq_u64(sum1, tVec);
                        uint64x2_t cmp2 = vcgtq_u64(sum2, tVec);

                        tmpVal1 = vsubq_u64(sum1, vandq_u64(cmp1, tVec));
                        tmpVal2 = vsubq_u64(sum2, vandq_u64(cmp2, tVec));

                        vst1q_u64(&scheme.context.polyTBootStrap[j][tN + n], tmpVal1);
                        vst1q_u64(&scheme.context.polyTBootStrap[j][tN + n + 2], tmpVal2);
                    }
                }
            }
        }
    }

    for (int i = 0; i < blockLen; i++) {
        scheme.context.TINTTAndEqual(scheme.context.polyTBootStrap[i], scheme.context.T);
        scheme.context.TAddAndEqual(scheme.context.polyTBootStrap[i], scheme.context.halfTPolyT, scheme.context.T);
    }
    // scheme.context.BaseConvertT2PQl(scheme.context.polyPQPool, scheme.context.polyTBootStrap, levelQ);
    // scheme.context.TSubAndEqual(scheme.context.polyPQPool, scheme.context.halfTPolyPQ, levelQ, k);
    // // return scheme.context.BaseConvertPQl2Ql(scheme.context.polyPQPool, levelQ);
    // return scheme.context.polyPQPool;

    scheme.context.BaseConvertT2PQl(result, scheme.context.polyTBootStrap, levelQ);
    scheme.context.TSubAndEqual(result, scheme.context.halfTPolyPQ, levelQ, k);
    // return scheme.context.BaseConvertPQl2Ql(scheme.context.polyPQPool, levelQ);
    return result;
}

/*
 * @brief ：执行加法操作
 *          该函数对两个多项式的相应部分进行加法操作，并更新第一个多项式 `a`。执行过程包括对密文的前 k 部分和后 l 部分分别进行加法操作。
 * 
 * @param a：第一个多项式，作为加法的结果。
 * @param b：第二个多项式，参与加法运算。
 * @param l：对应模数Q长度。
 * @param k：对应模数P长度。
 */
void EncodingMatrix::AddAndEqual(uint64_t *a, uint64_t *b, long l, long k)
{
    int logN = scheme.context.logN;
    for (long i = 0; i < k; i++) {
        uint64_t *ai = a + (i << logN);
        uint64_t *bi = b + (i << logN);
        scheme.context.PiAddAndEqual(ai, bi, i);
    }

    for (long i = k; i < l + k; i++) {
        uint64_t *ai = a + (i << logN);
        uint64_t *bi = b + (i << logN);
        scheme.context.QiAddAndEqual(ai, bi, i - k);
    }
}

/*
 * @brief ：执行 NTT 变换
 *          该函数对输入的多项式进行 NTT 变换，并更新多项式 `a`。对于前 k 部分和后 l 部分分别执行不同的变换操作。
 * 
 * @param a：输入的多项式。
 * @param l：对应模数Q长度。
 * @param k：对应模数P长度。
 */
void EncodingMatrix::NTTAndEqual(uint64_t *a, long l, long k)
{
    int logN = scheme.context.logN;
    for (long index = 0; index < k; index++) {
        uint64_t *ai = a + (index << logN);
        scheme.context.PiNTTAndEqual(ai, index);
    }

    for (long index = k; index < l + k; index++) {
        uint64_t *ai = a + (index << logN);
        scheme.context.QiNTTAndEqual(ai, index - k);
    }
}

/*
 * @brief ：根据给定的索引进行多项式排列
 *          该函数根据给定的索引对输入的多项式进行排列，生成输出多项式。
 * 
 * @param polOut：输出多项式。
 * @param polIn：输入多项式。
 * @param index：排列索引数组。
 * @param l：对应模数Q长度。
 * @param k：对应模数P长度。
 */
void EncodingMatrix::PermuteNTTWithIndex(uint64_t *polOut, uint64_t *polIn, uint64_t *index, int l, int k)
{
    int logN = scheme.context.logN;
    for (long i = 0; i < l + k; i++) {
        uint64_t *polIni = polIn + (i << logN);
        uint64_t *polOuti = polOut + (i << logN);
        scheme.context.IPermuteNTTWithIndex(polOuti, polIni, index);
    }
}

/*
 * @brief ：执行同态旋转操作的扩展
 *          该函数首先获取旋转密钥，然后执行内积计算，进行 NTT 变换，最后进行旋转并将结果返回。该操作用于同态加密中的加密旋转。
 * 
 * @param result：存储旋转操作结果的输出密文。
 * @param cipher：输入密文对象。
 * @param digits：多项式的旋转数字。
 * @param index：旋转操作的索引。
 * @param addFirst：是否在旋转前先加法。
 */
void EncodingMatrix::EvalFastRotationExt(Ciphertext &result, Ciphertext &cipher, vector<uint64_t *> digits,
                                         int32_t index, bool addFirst)
{
    uint64_t N = scheme.context.N;
    int K = scheme.context.K;
    int logN = scheme.context.logN;
    scheme.AddLeftRotKey(secretKey, index);
    SwitchKey switchKey = scheme.leftRotKeyMap.at(index);
    uint64_t *bx = new uint64_t[(cipher.l + K) * N]();
    if (addFirst) {
        std::copy(cipher.bx, cipher.bx + cipher.l * N, bx + K * N);
    } else {
        std::copy(cipher.bx, cipher.bx + (cipher.l + K) * N, bx);
    }
    result.ax = new uint64_t[(cipher.l + K) * N]();
    result.bx = new uint64_t[(cipher.l + K) * N]();
    uint64_t *axRotate = new uint64_t[(cipher.l + K) << logN]();
    uint64_t *bxRotate = new uint64_t[(cipher.l + K) << logN]();
    axRotate = InnerProduct(cipher.l, digits, switchKey.rax);
    bxRotate = InnerProduct(cipher.l, digits, switchKey.rbx);
    NTTAndEqual(axRotate, cipher.l, K);
    NTTAndEqual(bxRotate, cipher.l, K);
    AddAndEqual(bxRotate, bx, cipher.l, K);
    uint64_t galElInv = scheme.context.FastPow(5, index, scheme.context.M);
    uint64_t *ntt_index = new uint64_t[N]();
    scheme.context.PermuteNTTIndex(galElInv, ntt_index);
    PermuteNTTWithIndex(result.ax, axRotate, ntt_index, cipher.l, K);
    PermuteNTTWithIndex(result.bx, bxRotate, ntt_index, cipher.l, K);
    delete[] axRotate;
    delete[] bxRotate;
    delete[] ntt_index;
    delete[] bx;
}

/*
 * @brief ：计算并返回旋转的有效索引
 *          该函数用于根据给定的索引和槽数，计算旋转后的有效索引。主要用于优化旋转操作的索引计算，适用于槽数为 2 的幂的情况。
 * 
 * @param index：旋转操作的索引。
 * @param slots：密文的槽数。
 * 
 * @return 返回旋转后的有效索引。
 */
uint32_t EncodingMatrix::ReduceRotation(int32_t index, uint32_t slots)
{
    int32_t islots = int32_t(slots);

    // if slots is a power of 2
    if ((slots & (slots - 1)) == 0) {
        int32_t n = std::log2(slots);
        if (index >= 0) {
            return index - ((index >> n) << n);
        }
        return index + islots + ((int32_t(std::fabs(index)) >> n) << n);
    }
    if (islots == 0) {
        throw std::runtime_error("Divisor mod by zero is not allowed.");
    }
    return (islots + index % islots) % islots;
}

/*
 * @brief ：对密文执行标量乘法
 *          该函数执行对密文 `a` 的每个元素进行乘大数P，并更新结果 `res`。
 * 
 * @param res：存储标量乘法结果的数组。
 * @param a：输入多项式。
 * @param l：对应模数Q长度。
 * @param k：对应模数P长度。
 */
void EncodingMatrix::MulScalarP(uint64_t *res, uint64_t *a, long l)
{
    int slice;
    for (int i = 0; i < l * scheme.context.N; i++) {
        slice = i / scheme.context.N;
        for (int k = 0; k < scheme.context.K; k++) {
            if (k == 0) {
                MulModBarrett(res[i], a[i], scheme.context.pVec[k], scheme.context.qVec[slice],
                              scheme.context.twoqVec[slice], scheme.context.qrVecHigh[slice], scheme.context.qrVecLow[slice]);
            } else {
                MulModBarrett(res[i], res[i], scheme.context.pVec[k], scheme.context.qVec[slice],
                              scheme.context.twoqVec[slice], scheme.context.qrVecHigh[slice], scheme.context.qrVecLow[slice]);
            }
        }
    }
}

/*
 * @brief ：执行加法扩展操作（就地修改）
 *          该函数对两个密文对象 `cipher` 和 `ext` 进行加法操作，并将结果就地更新到 `cipher` 中。操作分为前 K 部分和后 l 部分分别进行加法计算。
 * 
 * @param cipher：输入密文，执行加法操作后将更新其值。
 * @param ext：扩展密文，将与 `cipher` 进行加法操作。
 */
void EncodingMatrix::EvalAddExtInPlace(Ciphertext &cipher, Ciphertext &ext)
{
    long K = scheme.context.K;
    long logN = scheme.context.logN;
    for (long i = 0; i < K; i++) {
        uint64_t *a1i = cipher.ax + (i << logN);
        uint64_t *b1i = cipher.bx + (i << logN);

        uint64_t *a2i = ext.ax + (i << logN);
        uint64_t *b2i = ext.bx + (i << logN);
        scheme.context.PiAddAndEqual(a1i, a2i, i);
        scheme.context.PiAddAndEqual(b1i, b2i, i);
    }

    for (long i = K; i < K + cipher.l; i++) {
        uint64_t *a1i = cipher.ax + (i << logN);
        uint64_t *b1i = cipher.bx + (i << logN);

        uint64_t *a2i = ext.ax + (i << logN);
        uint64_t *b2i = ext.bx + (i << logN);
        scheme.context.QiAddAndEqual(a1i, a2i, i - K);
        scheme.context.QiAddAndEqual(b1i, b2i, i - K);
    }
}

/*
 * @brief ：执行 PQ 重缩放并解密操作
 *          该函数对密文进行 PQ 重缩放，并将结果解密。
 * 
 * @param secretKeyPQ：私钥，解密所需的密钥。
 * @param ct：输入密文，进行解密操作。
 * @param s：输出标签，输出解密结果时的描述。
 */
void EncodingMatrix::PQRescaleAndDecrypt(SecretKey &secretKeyPQ, Ciphertext &ct, string s)
{
    long K = scheme.context.K;
    uint64_t *ctax = new uint64_t[(ct.l + K) * ct.N]();
    uint64_t *ctbx = new uint64_t[(ct.l + K) * ct.N]();
    std::copy(ct.ax, ct.ax + (ct.l + K) * ct.N, ctax);
    std::copy(ct.bx, ct.bx + (ct.l + K) * ct.N, ctbx);
    for (long index = 0; index < scheme.context.K; index++) {
        uint64_t *ai = ctax + (index << scheme.context.logN);
        uint64_t *bi = ctbx + (index << scheme.context.logN);
        scheme.context.PiINTTAndEqual(ai, index);
        scheme.context.PiINTTAndEqual(bi, index);
    }

    for (long index = K; index < ct.l + K; index++) {
        uint64_t *ai = ctax + (index << scheme.context.logN);
        uint64_t *bi = ctbx + (index << scheme.context.logN);
        scheme.context.QiINTTAndEqual(ai, index - K);
        scheme.context.QiINTTAndEqual(bi, index - K);
    }

    uint64_t *ax = new uint64_t[ct.l * ct.N]();
    uint64_t *bx = new uint64_t[ct.l * ct.N]();
    // ax = scheme.context.BaseConvertPQl2Ql(ctax, ct.l);
    // bx = scheme.context.BaseConvertPQl2Ql(ctbx, ct.l);
    scheme.context.BaseConvertPQl2Ql(ax, ctax, ct.l);
    scheme.context.BaseConvertPQl2Ql(bx, ctbx, ct.l);
    scheme.context.NTTAndEqual(ax, ct.l);
    scheme.context.NTTAndEqual(bx, ct.l);
    Ciphertext result = Ciphertext(ax, bx, ct.N, ct.slots, ct.l, ct.scale);
    scheme.ReScaleByAndEqual(result, 1);
    complex<double> *vec = scheme.Decrypt(secretKeyPQ, result);
    cout << s << ":" << endl;
    for (int i = 0; i < 5; i++) {
        cout << vec[i + 1024] << endl;
    }
    cout << endl;
    // delete [] ax;
    // delete [] bx_CAE;
    // delete []ctax;
    // delete []ctbx;
}

/*
 * @brief ：执行 PQ 解密操作。该函数对密文进行 PQ 解密。
 * 
 * @param secretKeyPQ：私钥，解密所需的密钥。
 * @param ct：输入密文，进行解密操作。
 * @param s：输出标签，输出解密结果时的描述。
 */
void EncodingMatrix::PQDecrypt(SecretKey &secretKeyPQ, Ciphertext &ct, string s)
{
    long K = scheme.context.K;
    uint64_t *ctax = new uint64_t[(ct.l + K) * ct.N]();
    uint64_t *ctbx = new uint64_t[(ct.l + K) * ct.N]();
    std::copy(ct.ax, ct.ax + (ct.l + K) * ct.N, ctax);
    std::copy(ct.bx, ct.bx + (ct.l + K) * ct.N, ctbx);
    for (long index = 0; index < scheme.context.K; index++) {
        uint64_t *ai = ctax + (index << scheme.context.logN);
        uint64_t *bi = ctbx + (index << scheme.context.logN);
        scheme.context.PiINTTAndEqual(ai, index);
        scheme.context.PiINTTAndEqual(bi, index);
    }

    for (long index = K; index < ct.l + K; index++) {
        uint64_t *ai = ctax + (index << scheme.context.logN);
        uint64_t *bi = ctbx + (index << scheme.context.logN);
        scheme.context.QiINTTAndEqual(ai, index - K);
        scheme.context.QiINTTAndEqual(bi, index - K);
    }

    uint64_t *ax = new uint64_t[ct.l * ct.N]();
    uint64_t *bx = new uint64_t[ct.l * ct.N]();
    // ax = scheme.context.BaseConvertPQl2Ql(ctax, ct.l);
    // bx = scheme.context.BaseConvertPQl2Ql(ctbx, ct.l);
    scheme.context.BaseConvertPQl2Ql(ax, ctax, ct.l);
    scheme.context.BaseConvertPQl2Ql(bx, ctbx, ct.l);
    scheme.context.NTTAndEqual(ax, ct.l);
    scheme.context.NTTAndEqual(bx, ct.l);
    Ciphertext result = Ciphertext(ax, bx, ct.N, ct.slots, ct.l, ct.scale);
    complex<double> *vec = scheme.Decrypt(secretKeyPQ, result);
    cout << s << ":" << endl;
    for (int i = 0; i < 5; i++) {
        cout << vec[i] << endl;
    }
    cout << endl;
    // delete []ax;
    // delete []bx_CAE;
    // delete []ctax;
    // delete []ctbx;
}

/*
 * @brief ：执行 Q 重缩放并解密操作
 
 * @param secretKeyQ：私钥，解密所需的密钥。
 * @param ct：输入密文，进行解密操作。
 * @param s：输出标签，输出解密结果时的描述。
 */
void EncodingMatrix::QRescaleAndDecrypt(SecretKey &secretKeyQ, Ciphertext &ct, string s)
{
    Ciphertext tmp = ct;
    scheme.ReScaleByAndEqual(tmp, 1);
    complex<double> *vec = scheme.Decrypt(secretKeyQ, tmp);
    cout << s << ":" << endl;
    for (int i = 0; i < 10; i++) {
        cout << vec[i] << " , ";
    }
    cout << endl;
}

/*
 * @brief ：执行 Q 解密操作
 * 
 * @param key：私钥，解密所需的密钥。
 * @param ct：输入密文，进行解密操作。
 * @param s：输出标签，输出解密结果时的描述。
 */
void EncodingMatrix::QDecrypt(SecretKey &key, Ciphertext &ct, string s)
{
    complex<double> *vec = scheme.Decrypt(key, ct);
    cout << s << ":" << endl;
    for (int i = 0; i < 20; i++) {
        cout << vec[i] << endl;
    }
    cout << endl;
}

/*
 * @brief ：执行同态乘法操作
 * 
 * @param result：存储乘法结果的密文对象。
 * @param cipher：输入密文，将与明文进行乘法运算。
 * @param pt：输入明文，用于与密文进行乘法运算。
 */
void EncodingMatrix::EvalMultExt(Ciphertext &result, Ciphertext &cipher, Plaintext pt)
{
    long L = cipher.l;
    long K = scheme.context.K;
    long logN = scheme.context.logN;
    long N = 1 << logN;

    result = cipher;

    result.ax = new uint64_t[(L + K) << logN]();
    result.bx = new uint64_t[(L + K) << logN]();

    for (long i = 0; i < K; i++) {
        uint64_t *ai = cipher.ax + (i << logN);
        uint64_t *bi = cipher.bx + (i << logN);
        uint64_t *mi = pt.mx + (i << logN);
        uint64_t *reai = result.ax + (i << logN);
        uint64_t *rebi = result.bx + (i << logN);

        scheme.context.PiMul(reai, ai, mi, i);
        scheme.context.PiMul(rebi, bi, mi, i);
    }
    for (long i = K; i < L + K; i++) {
        uint64_t *ai = cipher.ax + (i << logN);
        uint64_t *bi = cipher.bx + (i << logN);
        uint64_t *mi = pt.mx + (i << logN);
        uint64_t *reai = result.ax + (i << logN);
        uint64_t *rebi = result.bx + (i << logN);
        scheme.context.QiMul(reai, ai, mi, i - K);
        scheme.context.QiMul(rebi, bi, mi, i - K);
    }

    result.scale = cipher.scale * cipher.scale;

    // return Ciphertext(axMid,bxMid, cipher.N, cipher.slots, cipher.l,cipher.scale * scheme.context.p);
}

/*
 * @brief ：对密文执行从系数到槽上的转化
 * 
 * @param A ：预生成的切换矩阵
 * @param cipher ：密文
 * @return Ciphertext ：返回转化后的结果密文
 */
Ciphertext EncodingMatrix::EvalCoeffsToSlots(const vector<vector<Plaintext>> &A, Ciphertext &cipher)
{
    int stop = -1;
    int flagRem = 0;
    int BASE_NUM_LEVELS_TO_DROP = 1;
    long N = scheme.context.N;
    long M = scheme.context.M;

    long K = scheme.context.K;
    long logN = scheme.context.logN;
    long slots = cipher.slots;

    int levelBudget = mParamsEnc[0];
    int layersCollapse = mParamsEnc[1];
    int remCollapse = mParamsEnc[2];
    int numRotations = mParamsEnc[3];
    int b = mParamsEnc[4];
    int g = mParamsEnc[5];
    int numRotationsRem = mParamsEnc[6];
    int bRem = mParamsEnc[7];
    int gRem = mParamsEnc[8];

    if (remCollapse != 0) {
        stop = 0;
        flagRem = 1;
    }

    std::vector<std::vector<int32_t>> rot_in(levelBudget);
    for (uint32_t i = 0; i < uint32_t(levelBudget); i++) {
        if (flagRem == 1 && i == 0) {
            rot_in[i] = std::vector<int32_t>(numRotationsRem + 1);
        } else {
            rot_in[i] = std::vector<int32_t>(numRotations + 1);
        }
    }

    std::vector<std::vector<int32_t>> rot_out(levelBudget);
    for (uint32_t i = 0; i < uint32_t(levelBudget); i++) {
        rot_out[i] = std::vector<int32_t>(b + bRem);
    }

    for (int32_t s = levelBudget - 1; s > stop; s--) {
        for (int32_t j = 0; j < g; j++) {
            rot_in[s][j] = ReduceRotation((j - int32_t((numRotations + 1) / 2) + 1) *
                                          (1 << ((s - flagRem) * layersCollapse + remCollapse)),
                                          slots);
        }

        for (int32_t i = 0; i < b; i++) {
            rot_out[s][i] = ReduceRotation((g * i) * (1 << ((s - flagRem) * layersCollapse + remCollapse)), M / 4);
        }
    }

    if (flagRem) {
        for (int32_t j = 0; j < gRem; j++) {
            rot_in[stop][j] = ReduceRotation((j - int32_t((numRotationsRem + 1) / 2) + 1), slots);
        }

        for (int32_t i = 0; i < bRem; i++) {
            rot_out[stop][i] = ReduceRotation((gRem * i), M / 4);
        }
    }

    Ciphertext result = cipher;
    for (int32_t s = levelBudget - 1; s > stop; s--) {
        if (s != levelBudget - 1) {
            scheme.ReScaleByAndEqual(result, BASE_NUM_LEVELS_TO_DROP);
        }
        auto digits = EvalFastRotationPrecompute(result);
        MulScalarP(result.ax, result.ax, result.l);
        MulScalarP(result.bx, result.bx, result.l);
        std::vector<Ciphertext> fastRotation(g);
        for (int32_t j = 0; j < g; j++) {
            if (rot_in[s][j] != 0) {
                fastRotation[j] = result;
                EvalFastRotationExt(fastRotation[j], result, digits, rot_in[s][j], true);
            } else {
                KeySwitchExt(fastRotation[j], result);
            }
        }
        Ciphertext outer;
        uint64_t *first = new uint64_t[result.l << logN]();
        for (int32_t i = 0; i < b; i++) {
            int32_t G = g * i;
            Ciphertext inner;
            EvalMultExt(inner, fastRotation[0], A[s][G]);
            for (int32_t j = 1; j < g; j++) {
                if ((G + j) != int32_t(numRotations)) {
                    Ciphertext tmp;
                    EvalMultExt(tmp, fastRotation[j], A[s][G + j]);
                    EvalAddExtInPlace(inner, tmp);
                }
            }
            if (i == 0) {
                outer = inner;
                outer.ax = new uint64_t[(result.l + K) << logN]();
                outer.bx = new uint64_t[(result.l + K) << logN]();
                std::copy(inner.ax, inner.ax + (result.l + K) * result.N, outer.ax);
                std::copy(inner.bx, inner.bx + (result.l + K) * result.N, outer.bx);
            } else {
                if (rot_out[s][i] != 0) {
                    Ciphertext tmp = KeySwitchDownFirst(inner);
                    auto innerDigits = EvalFastRotationPrecompute(tmp);
                    EvalFastRotationExt(tmp, inner, innerDigits, rot_out[s][i], false);
                    EvalAddExtInPlace(outer, tmp);
                } else {
                    EvalAddExtInPlace(outer, inner);
                }
            }
        }
        result = KeySwitchDown(outer);
    }

    vector<Ciphertext> fastRotationRem(gRem);
    if (flagRem) {
        scheme.ReScaleByAndEqual(result, BASE_NUM_LEVELS_TO_DROP);
        auto digits = EvalFastRotationPrecompute(result);
        MulScalarP(result.ax, result.ax, result.l);
        MulScalarP(result.bx, result.bx, result.l);
        for (int32_t j = 0; j < gRem; j++) {
            if (rot_in[stop][j] != 0) {
                fastRotationRem[j] = result;
                EvalFastRotationExt(fastRotationRem[j], result, digits, rot_in[stop][j], true);
            } else {
                KeySwitchExt(fastRotationRem[j], result);
            }
        }
        Ciphertext outer;
        for (int32_t i = 0; i < bRem; i++) {
            Ciphertext inner;
            int32_t GRem = gRem * i;
            EvalMultExt(inner, fastRotationRem[0], A[stop][GRem]);
            for (int32_t j = 1; j < gRem; j++) {
                if ((GRem + j) != int32_t(numRotationsRem)) {
                    Ciphertext tmp;
                    EvalMultExt(tmp, fastRotationRem[j], A[stop][GRem + j]);
                    EvalAddExtInPlace(inner, tmp);
                }
            }

            if (i == 0) {
                outer = inner;
                outer.ax = new uint64_t[(result.l + K) << logN]();
                outer.bx = new uint64_t[(result.l + K) << logN]();
                std::copy(inner.ax, inner.ax + (result.l + K) * result.N, outer.ax);
                std::copy(inner.bx, inner.bx + (result.l + K) * result.N, outer.bx);
            } else {
                if (rot_out[stop][i] != 0) {
                    Ciphertext tmp = KeySwitchDownFirst(inner);

                    auto innerDigits = EvalFastRotationPrecompute(tmp);
                    EvalFastRotationExt(tmp, inner, innerDigits, rot_out[stop][i], false);
                    EvalAddExtInPlace(outer, tmp);
                } else {
                    EvalAddExtInPlace(outer, inner);
                }
            }
        }

        result = KeySwitchDown(outer);
    }
    scheme.ReScaleByAndEqual(result, 1);
    return result;
}

/*
 * @brief ：对密文执行从槽到系数上的转化
 * 
 * @param A ：预生成的切换矩阵
 * @param cipher ：密文
 * @return Ciphertext ：返回转化后的结果密文
 */
Ciphertext EncodingMatrix::EvalSlotsToCoeffs(const vector<vector<Plaintext>> &A, Ciphertext &cipher)
{
    int BASE_NUM_LEVELS_TO_DROP = 1;
    int levelBudget = mParamsDec[0];
    int layersCollapse = mParamsDec[1];
    int remCollapse = mParamsDec[2];
    int numRotations = mParamsDec[3];
    int b = mParamsDec[4];
    int g = mParamsDec[5];
    int numRotationsRem = mParamsDec[6];
    int bRem = mParamsDec[7];
    int gRem = mParamsDec[8];
    int flagRem = 0;
    long N = scheme.context.N;
    long logN = scheme.context.logN;
    long M = scheme.context.M;
    long K = scheme.context.K;
    long slots = N / 2;
    if (remCollapse != 0) {
        flagRem = 1;
    }

    vector<vector<int32_t>> rot_in(levelBudget);
    for (uint32_t i = 0; i < uint32_t(levelBudget); i++) {
        if (flagRem == 1 && i == uint32_t(levelBudget - 1)) {
            rot_in[i] = vector<int32_t>(numRotationsRem + 1);
        } else {
            rot_in[i] = vector<int32_t>(numRotations + 1);
        }
    }

    vector<vector<int32_t>> rot_out(levelBudget);
    for (uint32_t i = 0; i < uint32_t(levelBudget); i++) {
        rot_out[i] = vector<int32_t>(b + bRem);
    }

    for (int32_t s = 0; s < levelBudget - flagRem; s++) {
        for (int32_t j = 0; j < g; j++) {
            rot_in[s][j] =
                    ReduceRotation((j - int32_t((numRotations + 1) / 2) + 1) * (1 << (s * layersCollapse)), M / 4);
        }

        for (int32_t i = 0; i < b; i++) {
            rot_out[s][i] = ReduceRotation((g * i) * (1 << (s * layersCollapse)), M / 4);
        }
    }

    if (flagRem) {
        int32_t s = levelBudget - flagRem;
        for (int32_t j = 0; j < gRem; j++) {
            rot_in[s][j] =
                    ReduceRotation((j - int32_t((numRotationsRem + 1) / 2) + 1) * (1 << (s * layersCollapse)), M / 4);
        }

        for (int32_t i = 0; i < bRem; i++) {
            rot_out[s][i] = ReduceRotation((gRem * i) * (1 << (s * layersCollapse)), M / 4);
        }
    }

    Ciphertext result = cipher;

    for (int32_t s = 0; s < levelBudget - flagRem; s++) {
        if (s != 0) {
            scheme.ReScaleByAndEqual(result, BASE_NUM_LEVELS_TO_DROP);
        }
        auto digits = EvalFastRotationPrecompute(result);
        MulScalarP(result.ax, result.ax, result.l);
        MulScalarP(result.bx, result.bx, result.l);
        std::vector<Ciphertext> fastRotation(g);
        for (int32_t j = 0; j < g; j++) {
            if (rot_in[s][j] != 0) {
                fastRotation[j] = result;
                EvalFastRotationExt(fastRotation[j], result, digits, rot_in[s][j], true);
            } else {
                KeySwitchExt(fastRotation[j], result);
            }
        }
        Ciphertext outer;
        uint64_t *first = new uint64_t[result.l << logN]();
        for (int32_t i = 0; i < b; i++) {
            int32_t G = g * i;
            Ciphertext inner;
            EvalMultExt(inner, fastRotation[0], A[s][G]);
            for (int32_t j = 1; j < g; j++) {
                if ((G + j) != int32_t(numRotations)) {
                    Ciphertext tmp;
                    EvalMultExt(tmp, fastRotation[j], A[s][G + j]);
                    EvalAddExtInPlace(inner, tmp);
                }
            }
            if (i == 0) {
                outer = inner;
                outer.ax = new uint64_t[(result.l + K) << logN]();
                outer.bx = new uint64_t[(result.l + K) << logN]();
                std::copy(inner.ax, inner.ax + (result.l + K) * result.N, outer.ax);
                std::copy(inner.bx, inner.bx + (result.l + K) * result.N, outer.bx);
            } else {
                if (rot_out[s][i] != 0) {
                    Ciphertext tmp = KeySwitchDownFirst(inner);
                    auto innerDigits = EvalFastRotationPrecompute(tmp);
                    EvalFastRotationExt(tmp, inner, innerDigits, rot_out[s][i], false);
                    EvalAddExtInPlace(outer, tmp);
                } else {
                    EvalAddExtInPlace(outer, inner);
                }
            }
        }
        result = KeySwitchDown(outer);
    }

    vector<Ciphertext> fastRotationRem(gRem);
    if (flagRem) {
        scheme.ReScaleByAndEqual(result, BASE_NUM_LEVELS_TO_DROP);
        auto digits = EvalFastRotationPrecompute(result);
        MulScalarP(result.ax, result.ax, result.l);
        MulScalarP(result.bx, result.bx, result.l);
        int32_t s = levelBudget - flagRem;
        for (int32_t j = 0; j < gRem; j++) {
            if (rot_in[s][j] != 0) {
                fastRotationRem[j] = result;
                EvalFastRotationExt(fastRotationRem[j], result, digits, rot_in[s][j], true);
            } else {
                KeySwitchExt(fastRotationRem[j], result);
            }
        }
        Ciphertext outer;
        for (int32_t i = 0; i < bRem; i++) {
            Ciphertext inner;
            int32_t GRem = gRem * i;
            EvalMultExt(inner, fastRotationRem[0], A[s][GRem]);
            for (int32_t j = 1; j < gRem; j++) {
                if ((GRem + j) != int32_t(numRotationsRem)) {
                    Ciphertext tmp;
                    EvalMultExt(tmp, fastRotationRem[j], A[s][GRem + j]);
                    EvalAddExtInPlace(inner, tmp);
                }
            }

            if (i == 0) {
                outer = inner;
                outer.ax = new uint64_t[(result.l + K) << logN]();
                outer.bx = new uint64_t[(result.l + K) << logN]();
                std::copy(inner.ax, inner.ax + (result.l + K) * result.N, outer.ax);
                std::copy(inner.bx, inner.bx + (result.l + K) * result.N, outer.bx);
            } else {
                if (rot_out[s][i] != 0) {
                    Ciphertext tmp = KeySwitchDownFirst(inner);
                    auto innerDigits = EvalFastRotationPrecompute(tmp);
                    EvalFastRotationExt(tmp, inner, innerDigits, rot_out[s][i], false);
                    EvalAddExtInPlace(outer, tmp);
                } else {
                    EvalAddExtInPlace(outer, inner);
                }
            }
        }

        result = KeySwitchDown(outer);
    }
    return result;
}

/*
 * @brief ：根据给定的槽数、预算和维度，计算并返回FFT参数
 * 
 * @param slots ：表示槽的数量
 * @param levelBudget ：表示级别预算
 * @param dim1 ：表示维度1的大小
 * @return std::vector<int32_t> ：返回包含压缩FFT参数的整型向量
 */
std::vector<int32_t> EncodingMatrix::GetCollapsedFFTParams(uint32_t slots, uint32_t levelBudget, uint32_t dim1)
{
    uint32_t logSlots = std::log2(slots);
    if (logSlots == 0) {
        logSlots = 1;
    }

    std::vector<uint32_t> dims = SelectLayers(logSlots, levelBudget);
    // Need to compute how many layers are collapsed in each of the level from the budget.
    // If there is no exact division between the maximum number of possible levels (log(slots)) and the
    // level budget, the last level will contain the remaining layers collapsed.
    int32_t layersCollapse = dims[0];
    int32_t remCollapse = dims[2];

    bool flagRem = remCollapse != 0;

    uint32_t numRotations = (1 << (layersCollapse + 1)) - 1;
    uint32_t numRotationsRem = (1 << (remCollapse + 1)) - 1;

    // Computing the baby-step b and the giant-step g for the collapsed layers for decoding.
    int32_t g;
    if (dim1 == 0 || dim1 > numRotations) {
        if (numRotations > 7) {
            g = (1 << (int32_t(layersCollapse / 2) + 2));
        } else {
            g = (1 << (int32_t(layersCollapse / 2) + 1));
        }
    } else {
        g = dim1;
    }
    int32_t b = (numRotations + 1) / g;

    int32_t bRem = 0;
    int32_t gRem = 0;
    if (flagRem) {
        if (numRotationsRem > 7) {
            gRem = (1 << (int32_t(remCollapse / 2) + 2));
        } else {
            gRem = (1 << (int32_t(remCollapse / 2) + 1));
        }
        bRem = (numRotationsRem + 1) / gRem;
    }

    // If this return statement changes then CKKS_BOOT_PARAMS should be altered as well
    return { int32_t(levelBudget),     layersCollapse, remCollapse, int32_t(numRotations), b, g,
             int32_t(numRotationsRem), bRem,           gRem };
}

/*
 * @brief ：对密文执行同态乘法运算
 * 
 * @param cipher ：输入密文
 * @param pt ：输入明文
 * @return Ciphertext ：返回经过扩展乘法后的密文
 */
Ciphertext EncodingMatrix::EvalMultExt(Ciphertext &cipher, Plaintext pt)
{
    long L = cipher.l;
    long K = scheme.context.K;
    long logN = scheme.context.logN;
    long N = 1 << logN;

    uint64_t *axMid = new uint64_t[(L + K) << logN]();
    uint64_t *bxMid = new uint64_t[(L + K) << logN]();

    for (long i = 0; i < K; i++) {
        uint64_t *ai = cipher.ax + (i << logN);
        uint64_t *bi = cipher.bx + (i << logN);
        uint64_t *mi = pt.mx + (i << logN);
        uint64_t *reai = axMid + (i << logN);
        uint64_t *rebi = bxMid + (i << logN);

        scheme.context.PiMul(reai, ai, mi, i);
        scheme.context.PiMul(rebi, bi, mi, i);
    }
    for (long i = K; i < L + K; i++) {
        uint64_t *ai = cipher.ax + (i << logN);
        uint64_t *bi = cipher.bx + (i << logN);
        uint64_t *mi = pt.mx + (i << logN);
        uint64_t *reai = axMid + (i << logN);
        uint64_t *rebi = bxMid + (i << logN);
        scheme.context.QiMul(reai, ai, mi, i - K);
        scheme.context.QiMul(rebi, bi, mi, i - K);
    }

    return Ciphertext(axMid, bxMid, cipher.N, cipher.slots, cipher.l, cipher.scale * cipher.scale);
}
